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ABSTRACT 


The  power  spectral  technique*  of  generalized  harmonic  analysis 
were  applied  to  determine  the  phase  and  amplitude  of  the  transfer 
function  between  the  vertical  component  of  turbulence  and  the 
resultant  pitch  rate  of  a  small  swept  wing  jet  fighter.  The  power 
spectrum  analysis  was  based  on  actual  flight  test  data.   In  addition 
the  amplitude  and  phase  of  the  aircraft's  transfer  function  are 
determined  theoretically  by  two  different  methods  and  compared  with 
the  results  of  the  power  spectrum  analysis.   It  was  established  that 
the  power  spectrum  method,  as  applied  to  random  turbulence,  can  be 
used  to  determine  an  aircraft's  response  to  this  random  input  on  a 
continuous  basis  and  that  power  spectrum  techniques  have  considerable 
potential  value  in  the  study  of  aircraft  gust  loading. 
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INTRODUCTION 

Atmospheric  turbulence  and  the  resulting  aircraft  disturbances 
from  stabilized  flight  due  to  this  turbulence  are  random  and  thus 
aperiodic  and  non-dissipative.  These  disturbances  can  only  be  de- 
fined on  a  non-precise  basis  with  certain  probability  distributions 
and  average  values  since  they  are  random  in  nature.  This  leads  to 
the  application  of  probability  and  statistical  concepts,  which,  in 
combination  with  classical  Fourier  techniques,  form  a  general 
approach  for  the  analysis  of  problems  of  this  type.  This  technique 
is  termed  the  power  spectrum  method. 

Power  spectral  methods  were  derived  from  economic  time  studies 
at  the  beginning  of  this  century.  The  initial  applications  were  also 
directed  toward  turbulence  studies  and  the  Brownian  motion  of  parti- 
cles. During  World  War  II  power  spectral  methods  were  used  to  analyze 
problems  in  servomechanisms  and  communications;  however,  the  applica- 
tion of  power  spectral  analysis  to  the  field  of  aeronautics,  and  in 
particular  to  the  study  of  the  dynamic  response  of  aircraft  to  tur- 
bulence, is  relatively  new  and  still  in  the  early  stages  of  development, 

Other  recent  studies  in  the  field  of  aeronautics  involving 
power  spectral  techniques  have  been  applied  to: 

a.  Wind  tunnel  turbulence 

b.  Gust  loads  on  aircraft 

c.  Acoustic  Noise 

d.  Buffeting  loads 

e.  Electronic  noise  in  guidance  and  gunnery  systems 

In  particular,  considering  the  gust  loading  and  dynamic  response 
of  aircraft,  statistical  inference  is  a  realistic  approach  since  the 
problem  is  essentially  one  of  a  random  process  phenomena.  Perhaps  a 
more  meaningful  presentation  of  the  effectiveness  of  the  power  spec- 
trum method  can  be  made  if  the  earlier  discrete  gust  design  methods  are 
mentioned  in  comparison. 
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The  earliest  discrete  gust  approach  involved  the  assumptions  of 
a  rigid  airplane,  a  ramp  type  gust,  and  a  specified  maximum  vertical 
velocity.   Subsequent  discrete  gust  design  methods  were  adjusted  for 
different  wing  chords  and  pitch  characteristics  between  different  air- 
craft by  expressing  data  in  terms  of  the  mass  ratio,  zd  .  The  dis- 
crete gust  concept  was  based  primarily  on  the  premise  that  the  rela- 
tive loads  for  single  isolated  gusts  are  a  measure  of  the  relative 
loads  in  a  sequence  of  gusts. 

However,  in  recent  years  a  more  general  approach  to  turbulence 
analysis  has  been  needed  in  order  to  account  for  the  greater  variety 
and  degree  of  aircraft  flexibility,  mass  distribution,  configuration 
and  dynamic  response  characteristics,  as  well  as  the  varied  operating 
speeds  and  altitudes.  Therefore,  the  power  spectral  techniques  of 
generalized  harmonic  analysis  were  developed  for  application  to  tur- 
bulence and  dynamic  response  studies. 

As  conditions  changed  and  the  variety  of  aircraft  and  operating 
altitudes/airspeeds  increased,  the  concept  of  the  discrete  gust  design 
method  was  broadened  to  include  flexible  bodies  as  well  as  rigid 
bodies.  Although  the  flexible  aircraft  was  mathematically  flown 
through  discrete  gusts  of  varying  wavelengths,  the  discrete  gust  con- 
cept did  not  adequately  account  for  the  continuous  nature  of  atmos- 
pheric  turbulence,  particularly  as  applied  to  lightly  damped  dynamic 
systems.  Latest  design  philosophy  embodies  both  the  discrete  gust 
concept  and  the  power  spectral  approach,  with  the  design  initially 
analyzed  and  checked  using  the  discrete  gust  method  and  then  examined 
in  greater  detail  using  the  power  spectrum  analysis.   Further  develop- 
ment of  the  power  spectral  method  may  result  in  its  use  in  the  future 
entirely  for  atmospheric  turbulence  studies. 

The  main  advantages  of  the  power  spectral  approach  to  turbu- 
lence analysis  are  outlined  in  Ref.  1: 

1.  It  provides  a  more  realistic  representation  of  the  -•  - 
continuous  nature  of  atmospheric  turbulence. 

2.  It  allows  airplane  configurations  and  response  character- 
istics to  be  taken  into  account  in  a  rational  manner. 

3.  It  allows  more  rational  consideration  of  design  and 
operational  variations  such  as  configuration  changes,  mission 
changes  and  airplane  degrees  of  freedom. 
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Most  of  the  present  day  applications  of  power  spectral  methods 
are  concerned  with  the  analysis  of  the  magnitude  and  distribution  of 
the  random  data;  for  example,  the  intensity  of  atmospheric  turbulence 
and  its  variation  with  frequency.   However  this  thesis  considers  a 
more  unique  application  of  power  spectral  techniques  by  establishing 
the  transfer  function  between  pitch  rate  and  turbulence  for  the 
F-100C  airplane  based  upon  the  analysis  of  a  random  input,  w  ,  the 
vertical  component  of  atmospheric  turbulence.   This  application  goes 
one  step  beyond  the  analysis  of  random  data  by  considering  its  effect 
upon  the  dynamic  response  of  an  airframe.  The  particular  transfer 
function.  —  ,  is  considered  to  provide  a  basic  indication  of  the 
F-100  airplane' 8  dynamic  response  characteristics  in  a  turbulent  en- 
vironment. To  measure  the  accuracy  of  the  results  obtained  through 
the  power  spectral  analysis,  a  comparison  between  —  ,  as  deter- 
mined  from  this  method,  is  made  with  the  values  of  JL  ,  as  deter- 
mined  by  solving  the  longitudinal  equations  of  dynamic  motion  for 
the  F-100C  airplane. 
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DISCUSSION  OF  THE  CHARACTERISTICS  AND  ANALYSIS  OF  RANDOM  DATA 

Any  recorded  data  which  describe  a  physical  process  can  be 
broadly  classified  into  either  one  of  two  categories:   deterministic  or 
non-deterministic.   Deterministic  data  are  information  which  describe  a 
process  which  can  be  represented  and  predicted  with  reasonable  accuracy 
by  exact  mathematical  relationships.   Some  common  examples  of  this  type 
of  process  are  the  electrical  potential  across  a  condenser  as  it  dis- 
charges through  a  resistor  or  the  motion  of  a  satellite  in  a  known 
orbit  around  the  earth.   Examples  of  non-deterministic  data  are  the 
variation  of  an  aircraft's  vertical  velocity  in  turbulence  and  the 
variation  in  the  magnitude  of  the  lift  developed  on  an  airfoil  encounter- 
ing buffeting  in  rough  air.   Due  to  the  inherent  randomness  of  the 
turbulence  this  non-deterministic  data  cannot  be  described,  nor  can  its 
behavior  be  predicted  by  explicit  mathematical  relationships.   In 
practical  terms,  the  decision  as  to  whether  or  not  physical  data  are 
classified  as  being  deterministic  or  non-deterministic  (random)  is 
based  upon  the  ability  to  reproduce  the  data  repeatedly  through  con- 
trolled experiments.   If  the  experiment  or  phenomena  which  produced 
the  data  can  be  repeated  with  the  same  results,  within  the  tolerances 
of  experimental  error,  then  the  data  obtained  are  classified  as  deter- 
ministic.  Conversely,  if  identical  results  cannot  be  obtained  upon 
repeating  an  experiment,  then  the  data  are  considered  to  be  non- 
deterministic  or  random  in  nature. 

In  considering  random  data,  therefore,  it  is  expected  that  an 
particular  observation  producing  this  data  will  effectively  represent 
only  one  of  many  possible  outcomes  which  might  have  occured  as  a  result 
of  this  experiment.   For  example,  assume  that  the  lift  developed  on  an 
airfoil  during  buffeting  in  flight  is  recorded  as  a  function  of  time.   A 
specific  lift  time  history  will  be  obtained  as  shown  in  Fig.  1. 
However,  if  a  second  airfoil  of  identical  shape  is  flown  through  a 
similar  type  of  disturbance  simultaneously,  it  is  likely  that  a  different 
lift  time  history  would  result.   Therefore,  the  lift  time  history  for 
any  one  airfoil  encountering  a  particular  finite  disturbance  is  only 
one  example  of  an  infinitely  large  number  of  time  histories  which  might 
have  occurred.   If  the  random  buffeting  is  observed  over  a  finite  time 
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Figure  I.   Sample  records  showing  the  lift  developed  on  an  airfoil 
in  turbulence. 
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interval  and  is  graphically  or  mathematically  represented  by  a  single 
time  history,  then  this  time  history  is  referred  to  as  a  sample  record 
by  Bendat  and  Piersol  (2).   The  collection  of  all  possible  sample 
records  of  lift  versus  time  which  the  random  buffeting  might  have  pro- 
duced is  called  a  random  process  or  stochastic  process  in  Ref.  2. 
As  a  result,  a  sample  record  of  observed  data  for  a  random  physical 
occurrence  can  be  considered  to  be  only  one  physical  outcome  of  a 
random  process. 

The  properties  of  a  physical  experiment  or  occurrence  which 
are  described  by  a  random  process  can  be  determined  from  the  average 
values  of  the  experimental  data  as  computed  over  the  range  of  sample 
records  obtained.   A  collection  of  sample  functions  which  comprise  a 
random  process  is  illustrated  in  Fig.  2.   By  taking  the  instantaneous 
values  of  each  of  the  sample  functions  at  time  t  ,  summing  the  values 
and  dividing  by  the  number  of  sample  functions,  the  mean  value  of  the 
random  process  for  a  particular  time  t  may  be  computed.   In  addition, 
a  joint  moment  or  correlation  between  the  values  of  the  random  pro- 
cess at  two  different  times  can  be  computed  by  taking  the  average  of 
the  product  of  the  instantaneous  values  at  two  different  times  t  and 
t     +    y  .      The  product  obtained  in  the  manner  described  above  is  termed 
the  auto-correlation  function.   The  mean  value,  //*  (^i)   »  an<*  tne  auto- 
correlation function/^  (ti^t/" y)  are  defined  as: 

*nd  /?„  (tt)  t,  +7)  =Jt«rO    ft  Z  Xa  (t)  *X,  (t,  +  7) 


(2) 


According  to  Ref.  2,  the  random  process  is  stationary  if 
^/£/x(£i)     and/^  (T/  /"/  +7)   do  not  vary  as  time  t.  is  varied.   In  a 
more  general  sense,  when  the  properties  of  a  random  process,  as  com- 
puted over  short  intervals,  do  not  vary  significantly  from  one  interval 
to  the  next,  then  the  process  is  considered  to  be  stationary. 
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xjt) 


Figure  2.   Ensemble  of  sample  functions  forming  a  random  process 
(Reproduced  from  Re f.  2). 
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Applying  the  concept  of  stationarity  to  atmospheric  turbulence, 
it  has  been  shown  in  Ref.  1  that  the  conditions  of  stationarity  and 
homogeneity  apply  to  large  scale  patterns  of  air  motion  within  regions 
of  one  hundred  miles  and  time  durations  of  one  hour.   In  cases  of  severe 
turbulence,  such  as  encountered  in  a  thunderstorm,  these  time  and  dis- 
tance estimates  would  not  be  valid.   However,  the  time  duration  for 
homogeneous  and  stationary  conditions  may  still  be  almost  five  minutes 
in  length.   It  is  usually  possible  to  determine  the  properties  of  a 
stationary  random  process  by  computing  time  averages  over  specified 
sample  functions  in  the  ensembles  (collection  of  sample  functions). 

Consider  the  k   sample  for  the  random  process  illustrated  in 
Fig.  2.  The  mean  value,  jC/x\-~)  ,   and  the  autocorrelation  function, 
/^  {C?  .*&.)      ,  for  the  k   sample  function,  according  to  Ref.  2, are: 

Random  processes  may  be  classified  as  shown  in  Fig.  3,  repro- 
duced from  Ref.  2,  as  being  either  stationary  or  non-stationary,  and 
stationary  random  processes  may  be  further  subdivided  as  either 
ergodic  or  non-ergodic.  The  random  process  is  ergodic  if  it  is  sta- 
tionary and ;£// (-&)  an6/?y(y^)do   not  vary  when  computed  over 
different  sample  functions.   Therefore,  all  properties  of  ergodic  ran- 
dom processes  can  be  determined  by  computing  over  a  single  sample 
function.   In  actual  physical  phenomena,  as  outlined  in  Ref.  2,  most 
random  data,  including  atmospheric  turbulence,  are  ergodic  so  that 
the  properties  of  the  data  can  be  computed  from  time  averages  of 
individual  sample  records. 

Two  types  of  statistical  functions  are  used  to  describe  the 
basic  properties  of  random  data:   (a)  the  autocorrelation  function 
and  (b)  the  power  spectral  density  function.   The  autocorrelation 
function  and  the  power  spectral  density  function  provide  an  analysis 
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Figure  3.   Classifications  of  random  data  (Reproduced  from  Ref.  2.) 
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of  the  amplitude  of  the  random  data  with  respect  to  time  and  frequency 
respectively.  These  functions  are  Fourier  transform  pairs,  and  there- 
fore provide  an  analysis  of  the  data  in  two  different  formats. 

The  autocorrelation  function  for  random  data  provides  an  indi- 
cation of  the  general  dependence  of  the  amplitude  of  the  data  at  a 
particular  time  on  the  value  of  the  data  at  a  different  time.  An 
autocorrelation  measurement  can  he  used  for  detecting  deterministic 
data  which  may  he  hidden  or  obscured  in  a  random  background. 

As  an  extension  of  the  autocorrelation  function,  the  power 
spectral  density  function  describes  the  frequency  composition  of  the 
data  in  terms  of  the  spectral  density  of  its  mean  square  value. 

In  Chapter  IV  an  explicit  and  more  detailed  presentation  is 
included  which  outlines  both  the  analytical  and  digital  computer  steps 
involved  in  the  power  spectrum  analysis  as  used  to  determine  the  pitch 
response  of  the  F-100C  airplane  to  atmospheric  turbulence.  This 
analysis  will  form  a  mathematical  extension  of  the  autocorrelation 
function,  through  Fourier  transformation,  to  yield  the  power  spectral 
density  function. 
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AIRCRAFT  DESCRIPTION  AMD  FLIGHT  TEST  PROCEDURES 

The  basic  measurement  used  to  determine  atmospheric  turbulence 
in  flight  is  the  angle  of  attack  of  the  aircraft  as  measured  by  flow 
vanes.  The  flow  angles  measured  are  corrected  for  the  effects  of 
airplane  motion,  Chilton  (3).   In  this  particular  analysis  as  applied 
to  the  power  spectrum  method,  a  continuous  time  history  of  the  vertical 
component  of  the  true  gust  velocity  is  to  be  initially  determined. 
Although  the  lateral  and  longitudinal  components  of  turbulence  could 
have  been  measured  as  well,  the  vertical  component  has  the  most  signi- 
ficant effect  upon  the  aircraft's  pitch  dynamics,  and  in  the  interest 
of  simplification  only  this  vertical  portion  of  the  turbulence  will  be 
considered. 

The  vertical  component  of  the  gust  velocity,  from  Ref.  1,  is 
given  by: 

iVa  =  Vocv  -  vq  /y&g  dt  t  Jx  e  (5) 

The  equation  above  is  based  upon  the  following  assumptions  as  outlined 
in  Ref.  1: 

1.  All  disturbances  are  small  enough  to  allow  the  use 
of  the  angle  in  radians  in  place  of  the  sine  of  the  angle. 

2.  Boom  flexures  are  negligible. 

3.  The  effects  of  variation  in  upwash  on  vane  indica- 
tions are  negligible  or  allowed  for  by  calibration. 

In  conducting  the  flight  analysis,  the  average  values  of  the 
measured  parameters  (denoted  by  a  bar)  are  included  to  establish  the 
flight  measurements  as  increments  from  the  mean  values  for  the  measured 
flight  interval.  Therefore,  the  time  history  of  the  vertical  compo- 
nent of  the  gust  velocity,  from  Ref.  1,  is: 

t 

u/<j  rz  v  (<*v-*v)-v(e-e)+f  (Q^-Q^dt  +^{4-5)     (6) 

A  reading  interval  of  0.1  seconds  was  established  as  being  an 
optimum  value  for  this  analysis  after  considering  the  aircraft  speed  at 
which  the  data  were  obtained  (300  knots  TAS)  and  the  frequencies 
evaluated. 
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All  turbulence  date  used  for  this  analysis  were  obtained  by 
the  author  during  the  month  of  June  1967  while  participating  in  an 
industrial  tour  at  the  Flight  Test  Branch  of  NASA's  Ames  Research 
Center  at  Moffett  Field,  California.  The  flight  test  data  used  for 
this  thesis  are  only  a  portion  of  the  total  data  obtained  by  NASA  as 
part  of  a  larger  program  to  study  the  optimum  dynamic  characteristics 
and  response  to  turbulence  of  a  variable  stability  F-100C  airplane. 

The  aircraft  used  for  this  turbulence  study  and  its  principal 
dimensions  are  shown  in  Fig.  4,  from  Ref.  4.   It  is  a  North  American 
F-100C  which  is  assigned  to  the  National  Aeronautics  and  Space  Adminis- 
tration and  designated  as  NASA  aircraft  number  703.  The  F-100  is  a 
single  place,  low  wing  aircraft,  powered  by  a  J -57  series  turbojet 
engine  with  afterburner.  The  wing  and  horizontal  stabilizer  are  swept 
back  at  an  angle  of  45  degrees.  Low  speed  flight  and  stall  character- 
istics are  improved  by  the  use  of  longitudinal  leading  edge  slats  on 
the  standard  F-100  airplane,  and  good  longitudinal  control  is  provided 
by  an  all-moving  "slab"  tail.  The  controls  are  completely  power 
actuated.  Longitudinal  stick  forces  are  provided  by  a  simple  spring 
bungee  whose  leverage  on  the  stick  is  constant  at  all  speed  ranges 
except  between  a  Mach  number  of  0.8  and  0.94  where  the  leverage  is 
altered  by  a  Mach  number  sensor.  The  F-100  airplane  has  a  400  cycle, 
115  volt  A.C.  electrical  system  and  a  separate  28  volt  D.C.  system. 

The  parameters  measured  in  flight  are  given  by  Equation  (5), 
and  the  type  and  location  of  the  flight  Instruments  are  shown  in 
Table  I.  A  continuous  record  of  aircraft  response  was  necessary  in 
order  to  provide  a  complete  and  effective  record  of  the  turbulence  en- 
countered, therefore  an  oscillograph,  instead  of  a  photopanel,  was 
required.  All  data  were  recorded  on  a  26  channel  photographically 
recording  oscillograph  and  developed  on  Ektachrome  Aero  film  for 
analysis.  As  can  be  seen  in  Fig.  5,  the  accelerometers  were  positioned 
as  close  to  the  aircraft's  center  of  gravity  as  possible.  The  flow 
vane  used  to  measure-  angle  of  attack  was  located  on  the  top  of  the 
vertical  stabilizer  as  shown.  A  standard  pitot  static  head  was 
utilized  for  airspeed  measurement.   It  is  located  ahead  of  the  aircraft 
on  the  tip  of  the  F-100' s  standard  nose  boom. 
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AIRPLANE  GENERAL  INFORMATION 
Airplane  Dimensions 


T.O.  1F-100C-2-5 


AIRPLANE  DIMENSIONS 


NUI  C  Airplane  principal  dimensiom  are  taken  with  the  landing 
gear  struts  and  tires  inflated  to  the  correct  pressures.  Weights 
given  are  approximate  and  are  not  to  be  used  for  weight  and  bal- 
ance computation  or  flight  restriction. 


FUSELAGE  REFERENCE  LINE 


AREAS 

WING  (TOTAL) - - 

AILERON  (TOTAL  AFT  OF    HINGE  LINE) 
WING  SLATS 


25  PERCENT  CHORD  LINE 


HORIZONTAL  STABILIZER   (TOTAL) 
VERTICAL  STABILIZER 


RUDDER   (AFT  OF  HINGE  LINE) 
SPEED  BRAKE 


385.2  SO  FT 
38.6  SQ  FT 
—  37.7  SQ  FT 
—99.0  SQ  FT 
—45.2  SQ  FT 
6  3  SQ  FT 


18.8  FT 


WEIGHT 

TAKE-OFF  GROSS  WEIGHT  (clean  airplane) 

MAXIMUM  TAKE-OFF  GROSS  WEIGHT     


WING 


TYPE 


AIRFOIL  SECTION 

ROOT- 
TIP 

INCIDENCE  AT  ROOT 

INCIDENCE  AT  TIP 

ASPECT  RATIO 


WING 
STABILIZER 


-NACA  64A007 
-NACA  64A007 
0  DEO 
0  DEG 
-3.86 


— J4.2  FT   |— 


f  ^^  \  6.4  FT 

(fl       * 

12.4  FT — 

38.8  FT 

F-100C-J-00-18D 


1-12 


Figure  4.   Airplane  Dimensions 
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PARAMETER 
MEASURED 


IIPE  OF  INSfRUHENT 


LOCATION 
OF  INSTRUMENT 


INSTRUMENT 
LIMITS 


Airspeed 
V 


Standard  F-100  pitot- 
static  instrument  with 
Bourns  pressure  cell 


Nose  boom 


0-500  kts 


Angle  of 
Attack 


1*00  cycle  \.C.  NASA 
instrument 


Top  of  vertical 
stabilizer, 
station  525 


*10° 


Pitch  Angle 


Mi/Jieapolis  Honeywell 
attitude  gyro  with 
potentiometer  pickoff 


Left  link 
compartment 
Station  152 


tuf 


Vertical 
Acceleration 

o, 


Statham  linear  strain 
gage  accelerometer. 
Model  CS-8-3SO 


Wheel  well  on 
aircraft ! s 
longitudinal 
center line 
Station  315 


£l*g«s 


Pitch  Rate 


Spring  restrained  rate 
gyro  with  a  Microsyn 
signal  generator 


Left  link 
compartment 
Station  ll*3 


±  90  /sec 


Longitudinal 
Stick  Position 


Potentiometer  with  i*00 
cycle  A.C.  output 


Stick  linkage 
on  centerline 


^19°stick 
deflection 


Horizontal 
Stabilizer 
Position 


Linear  potentiometer 
pickoff  using  Mar kite 
infinite  resolution  oot. 


Station  1*75 


^5°,  -5° 

stabilizer 

position 


Altitude 
h 


Statham  pressure  cell 
with  strain  gage  pickoff 


Nose  compartment 
Station  55 


0-50,000  ft 


TABLE  I:   AIRCRAFT  INSTRUMENTATION:   NASA  AIRCRAFT  NUMBER  703 
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This  particular  F-100C  had  been  modified  by  NASA  for  use  as  a 
variable  stability  airplane  in  order  to  study  the  effects  of  the 
interaction  between  the  turbulence  and  the  longitudinal  stability 

characteristics.   A  planned  variation  in  C    and  C   was  developed 

"V*       mq 
by  NASA  engineers,  and  these  two  parameters  were    M  changed  to  pro- 
vide varying  degrees  of  pitch  response. 

From  basic  theory,  the  values  of  C    primarily  determine  the 

m^c 

natural  frequency  of  the  short  period  mode  and  are  a  major  factor  in 
determining  the  response  of  the  airframe  to  elevator  motions  and  gusts. 

The  stability  derivative,  C   ,  is  usually  referred  to  as  the 

m 

pitch  damping  derivative  and  is    "  critical  in  longitudinal  dynamics 

because  it  establishes  a  major  portion  of  the  damping  of  the  short 

period  mode.   This  damping  effect  comes  mostly  from  the  horizontal 

stabilizer.   The  appropriate  C    and  C   gain  settings  were  obtained 

m06       mq 
by  means  of  an  analog  computer  matching  M  technique  to  produce  the 

desired  combinations  of  the  short  period  frequency  and  damping  ratio. 

For  this  analysis,  by  varying  C    and  C   ,  a  resultant  value 

mc*      mq 

of  M^   equal  to  -14.47  and  M  equal  to  .3068   M  was  programmed  into 

q 

the  aircraft  to  produce  a  short  period  frequency  of  .6  cycles  per 
second  and  a  damping  ratio  equal  to  0.10. 

The  simplified  schematics  in  Figs.  5  and  6  show  how  the  variable 
stability  mode  of  this  NASA  F-100C  is  obtained.   Three  modes  of  control 
may  be  manually  selected  by  the  pilot:   a.   Basic  Mode,  b.   "Fly  By 
Wire"  Mode  and  c.   Variable  Stability  Mode  (includes  "Fly  By  Wire" 
mode) . 

In  the  variable  stability  mode  as  shown  in  Fig.  5,  a  signal 
from  both  the  angle  of  attack  vane  and  the  pitch  rate  gyro  are  received 
at  the  summing  amplifier.   The  output  from  the  summing  amplifier  is 
summed  with  the  pilot's  input  and  fed  into  a  servo  which  controls  the 
horizontal  stabilizer's  hydraulic  actuating  valve  to  provide  the 
desired  response  characteristics. 

A  detailed  schematic  of  the  pilot's  input  as  a  closed  loop 
system  is  shown  in  Fig.  6.   A  potentiometer  is  used  to  detect  control 
stick  inputs  with  a  sensitivity  such  that  the  basic  aircraft  stick-to- 
stabilizer  gearing  is  preserved.   A  torque  tube  between  the  stick  and 
the  hydraulic  valve  is  normally   mated  firmly  in  the  Basic  Aircraft  Mode, 
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In  the  "Fly  By  Wire"  or  Variable  Stability  Modes,  the  torque  tube  is 
mechanically  adjusted  by  the  system  to  have  a  greater  increment  of 
play.  This  adjustment  prevents  any  feedback  from  the  stabilizer  to 
the  stick;  and  allows  additional  commands  to  the  stabilizer  in  the 
Variable  Stability  mode. 

C   and  C   are  altered  by  the  variable  stability  system  in 
accordance  with  the  following  equations: 


s£-#yc>_  (7) 

C^  e,r  -    Cn^  6as/c  +£&*ji   »Ae~  AGn^  G,fJg) 

£""$  eff  =  C™9.  6as<c  ^\^£f^j(±fjj  (8) 

In  each  example,  the  basic  values  of  the  F-100' s  stability 
derivatives  are  modified  by  the  second  term  in  the  equation  whose 
value  is  altered  by  the  servo  system. 

Since  the  longitudinal  response  of  the  F-100  airplane  to  tur- 
bulence was  of  prime  importance,  all  flight  tests  were  conducted  with 
no  longitudinal  control  inputs,  however  the  wings  were  kept  level, 
while  data  were  being  taken,  by  applying  minor  lateral  stick  inputs. 
As  a  normal  procedure,  the  turbulence  was  entered  from  a  stabilized 
flight  condition  at  an  airspeed  of  300  knots  true,  which  is  a  realistic 
cruise  airspeed  for  the  F-100  airplane  at  the  test  altitudes  selected 
(13,000  to  15,000  feet).  Although  no  longitudinal  control  inputs 
were  desired,  both  stabilizer  position  and  control  stick  deflection 
were  checked  continually  when  analyzing  the  flight  test  data  to  insure 
the  validity  of  the  measured  aircraft  response. 

Most  of  the  turbulence  encountered  was  of  a  mountain  wave  or 
clear  air  turbulence  type.  The  leeward  side  of  Mount  Whitney  and  the 
neighboring  peaks  in  the  Sierra  Nevada  range  were  used  as  a  flight 
test  area. 
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BASIC  COMPUTATION  PROCEDURE 

As  in  any  statistical  approach  to  analyzing  data  and  predicting 
trends,  the  more  information  obtained,  the  more  valid  the  study.  This 
concept  was  tempered  by  the  selection  of  that  flight  test  data  which 
were  obtained  with  the  F-100C  programmed  for  the  most  realistic  short 
period  frequency.  Accordingly,  after  screening  the  available  test  data, 
NASA  Flight  26,  Run  8  was  selected  as  being  a  valid  and  optimum  flight 
record  for  the  power  spectrum  analysis.  On  this  flight,  the  F-100C 
was  programmed  for  a  short  period  frequency  of  .6  cycles  per  second. 
It  is  reiterated  that  the  intent  of  this  thesis  is  to  demonstrate  the 
application  of  the  power  spectrum  analysis  to  determine  the  dynamic 
response  of  an  aircraft  to  turbulence  rather  than  to  establish  any  new 
information  on  the  F-100C  airplane. 

The  data  for  this  flight  were  read  every  0.1  seconds,  and  the 
basic  parameters  as  outlined  in  Equation  (5),  were  recorded  on  stan- 
dard IBM  eighty  column  digital  computer  cards.   Due  to  the  large  quanti- 
ties of  data  involved  (six  inflight  parameters  read  at  920  Intervals 
each)  and  the  numerous  computations  required  by  the  power  spectral 
method,  a  digital  computer  is  necessary  to  effectively  process  the 
data  and  perform  the  analysis  efficiently.  The  average  values  of  the 
measured  flight  parameters,  denoted  by  a  bar  in  Equation  (6),  were 
initially  determined  by  means  of  a  simple  digital  operation,  and  w 
itself  was  computed.  The  power  spectral  analysis  of  the  measured 
flight  data  was  performed  by  a  digital  program  utilizing  an  IBM 
7094/7040  direct  coupled  system  computer  for  the  calculations  and 
graphing. 

Both  the  autocorrelation  function  and  the  power  spectral  den- 

sity  of  the  input,  w  (referred  to  as  series  2),  and  the  output,  O 

(referred  to  as  series  1) ,  were  determined  and  graphed  by  the  computer 

program.  The  desired  transfer  function,  —   was  then  computed  in 

9 
both  phase  and  amplitude  as  a  function  of  gust  frequency.  The  transfer 

function  was  based  upon  a  cross  correlation  and  cross  spectral  density 

* 
function  which  combined  the  average  values  of  Q        and  w  over  a 

© 

common  frequency  interval.  From  the  computer  analysis,  the  phase  and 
amplitude  of   -^-    were  graphed  on  a  cross  spectral  density  plot 

"9 
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versus  frequency.  To  check  the  accuracy  of  the  power  spectrum  method, 

Q 

the  theoretical  solution  for  the  phase  and  amplitude  of  — —•   was 
obtained  by  solving  the  longitudinal  dynamic  equations  of  motion  and 
the  pitching  moment  equation  about  the  lateral  axis.  The  phase  and 
amplitude  of  the  transfer  function  obtained  theoretically  were  also 
plotted  versus  frequency,  and  the  resulting  curves  from  both  methods 
were  compared. 
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DISCUSSION  OF  THE  POWER  SPECTRUM  METHOD 

A  detailed  analysis  of  the  steps  involved  in  the  digital  com- 
puter program  is  presented  in  Appendix  A;  however,  a  proper  apprecia- 
tion of  the  digital  computation  can  only  be  realized  if  the  mathe- 
matical theory  behind  the  power  spectrum  method  is  outlined  in  detail. 

As  discussed  earlier,  the  vertical  component  of  the  atmospheric 

turbulence,  w  ,  was  established  as  a  stationary,  ergodlc  form  of 

g 
random  data.   An  additional  stipulation  is  that  the  turbulence  should 

represent  a  Gaussian  process.   For  the  turbulence  to  be  Gaussian,  the 
fluctuations  of  the  turbulence  must  have  a  normal  probability  distri- 
bution.  References  1  and  5  stipulate  that  if  one  component  of  the 
turbulence,  w  ,  for  instance,  has  a  normal  distribution,  then  normality 
of  the  turbulence  can  be  assumed.   Figure  7,  which  is  reproduced  from 
Ref.  1,  shows  the  cumulative  frequency  distribution  of  w  as  determined 
by  NASA  from  a  thunderstorm  traverse.   The  plot  is  scaled  in  such  a 

manner  that  a  normal  distribution  of  w  will  result  in  a  straight  line. 

g 
The  data  were  obtained  from  a  200  second  record,  read  at  0.1  second 

intervals.  A  fitted  normal  distribution  is  shown  in  comparison  to 

the  data  points  plotted,  and  it  may  be  noted  that  there  is  a  good 

correspondence  between  the  fitted  normal  distribution  and  the  traverse 

data  points,  except  at  large  values  of  w  ,  therefore,  from  Ref.  1, 

g 
the  assumption  that  there  is  a  Gaussian  distribution  of  turbulence 

appears  to  be  reasonable.   It  is  reiterated  that  the  properties  of 
stationarity,  homogeneity,  and  isotropy,  in  addition  to  the  Gaussian 
assumptions,  are  implied.  These  properties  specify,  then,  an  invariance 
in  the  statistical  characteristics  of  turbulence  with  respect  to  the 
position  along  the  time  history  of  the  sampling  interval  and  the 
direction  through  the  turbulence. 

From  Ref.  6,  Taylor's  hypothesis  that  the  variation  in  gust 
intensity  existing  along  the  flight  path  at  any  given  instant  will 
remain  substantially  the  same  until  the  aircraft  has  traversed  the 
given  body  of  air,  is  utilized.  The  relative  length  of  the  body  of 
air  along  the  flight  path,  therefore  is  large  compared  with  the  dis- 
tance traversed  in  the  reaction  time  of  the  aircraft. 
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Figure  7.  -  Probability  of  equaling  or  exceeding 

given  values  of  vertical  gust  velocity 
for  thunderstorm  (traverse  A). 
(Reproduced  from  Ref .  1) 
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With  the  mathematical  assumptions  established,  a  detailed 
analysis  of  the  mathematical  steps  behind  the  power  spectrum  analysis 
can  now  be  presented  in  conjunction  with  appropriate  sequential 
reference  to  the  computer  program  which  is  detailed  in  Appendix  A. 

A  block  diagram  (Fig.  17)  showing  the  steps  involved  in  the 
computer  analysis  and  the  appropriate  mathematical  expression  for 
each  step  is  presented  in  Appendix  A  in  order  to  provide  a  complete 
picture  of  the  entire  digital  program  involved  in  the  power  spectrum 
analysis. 

0 

The  magnitude  of  S     and  w  ,  as  computed  from  Equation  (6) , 
is  tabulated  at  0.1  second  intervals  on  the  computer  printout  and 
referred  to  as  series  1  and  series  2  respectively.  The  values  of  the 
data  are  read  from  left  to  right  in  each  row  in  chronological  order. 
The  autocorrelation  function  of  0     and  the  autocorrelation  function 
of  w  are  computed  separately  by  taking  the  product  of  0     or  w  at 
time  t  and  at  time  t  ♦  7,  and  averaging  these  values  over  a  total 
observation  time,  T.  The  resulting  average  product  will  approach  an 
exact aatfoeerrelation  function  (also  referred  to  as  the  autocovariance 
function)  as  T  approaches  infinity. 

In  equation  form  the  autocorrelation  function  is  expressed  as: 


7" 


(9) 

7-—*-  <»   Qy 

where  X (£)     represents  the  value  of  either  O    or  w  at  time  t,  and 

X(t+  y)    is  the  value  of  6  or  w  after  a  time  shift  whose  magnitude 

S 
is  determined  by  the  value  of  y   . 

The  autocorrelation  functions  of  Q     and  w  are  each  computed 

g 
for  progressive  values  of  y  ,  that  vary  in  0.1  second  increments, 

from  0.1  seconds  to  9.99  seconds.  The  autocorrelation  function  for 

S        is  tabulated  on  page  77  of  Appendix  A  and  plotted  versus  time 

on  page   78   also.   Similarly  /?*  for  w  is  tabulated  on  page  82 

and  plotted  versus  time  on  page  83  in  Appendix  A. 

The  quantity  rfyfe)  ,   as  defined  above,  is  a  real  valued  function 

with  a  maximum  value  at  y   -  0.  Also /?x(-7^RJy),   and  Rx(o)>  //?x(y)/ 

for  all  y  .  The  mean  value  of  x(t)  can  be  determined  from  Iff  (o£\   , 

i.e.,  it  is  equal  to  the  positive  square  root  of  the  autocorrelation 
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function,  as  the  time  displacement  increases  toward  infinity  (Ref.  2). 
In  addition,  the  mean  square  value  is  equal  to  the  autocorrelation 
function  at  zero  time  displacement,  (  jf   -  0) ;  i.e.,  &    —  ^x  ^°' 

As  an  example  of  the  usefulness  of  the  autocorrelation  function, 
a  simple  autocorrelogram  (plot  of  an  autocorrelation  function)  is 
illustrated  in  Figure  8. 


*M 


Figure  8.  -  Autocorrelogram  for  a  Sine  Wave. 
(Reproduced  from  Ref.  2) 


Although  phase  information  is  not  available  with  this  plot,  it 
can  be  noted  that  the  autocorrelation  function  for  the  sine  wave  will 
persist  over  all  time  displacements  since  it  represents  deterministic 
data.  Conversely,  an  autocorrelogram  for  narrow  band  random  noise, 
(non-deterministic  data),  will  diminish  to  zero  for  large  time  dis- 
placements, assuming  a  mean  value  of  zero.   This  is  illustrated  in 
Fig.  9 


Figure  9.  -  Autocorrelogram  for  Random  Noise 
(Reproduced  from  Ref.  2) 


A  combined  plot  of  the  autocorrelation  function  of  the  sine 
wave  superimposed  with  the  autocorrelation  function  of  the  random 
noise  will  appear  as  in  Figure  10. 
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Figure  10.  -  Combined  Autocorrelograms  for 
a  Sine  Wave  and  Random  Noise 
(Reproduced  from  Ref .  2) 


The  deterministic  data  (sine  wave)  can  be  detected  from  the 
random  background  of  the  wide  band  noise.  This  is  the  principle 
value  of  the  autocorrelation  function  in  that  the  random  and  non- 
random  data  may  be  separated  into  their  respective  categories  through 
analysis  of  their  combined  autocorrelograms. 

o 

The  autocorrelation  function  of  B   ,  series  1,  is  noted  to 
diminish  with  increasing  values  of  time,  due  to  the  inherent  longi- 
tudinal damping  characteristics  of  the  F-100  airplane.  This  is  a 
typical  response  of  a  damped  second  order  system  to  a  random  input. 
Refer  to  pages  77  and  78  of  Appendix  A. 

It  is  also  to  be  noted  that  the  autocorrelation  function  of 

w  ,  series  2,  is  also  a  maximum  at  J  ■  0  and  decreases  to  a  minimum 

g 
value  for  increasing  values  of  y   (pages   82   and   83  ,  Appendix  A) . 

This  is  the  expected  characteristic  of  the  autocorrelation  function 

for  random  data. 

After  establishing  each  autocorrelation  function,  /?g   and 

A^ft  ,  the  power  spectral  estimates  of  both  £  and  w  are  computed 

<f  g 

separately  and  tabulated  on  pages  79  and   84  of  Appendix  A.  The 

power  spectral  estimate  is  simply  a  Fourier  transform  of  the  auto- 
correlation function  and  describes  the  general  frequency  composition 
of  the  data  in  terms  of  the  spectral  density  of  its  mean  square  value, 

(//  (Refs.  2  and  7).  The  time  mean  square  value  provides  a 

measure  of  the  disturbance  energy  per  unit  time  and  thus  has  been 
referred  to  as  power  in  the  study  of  turbulence  and  has  been  termed 
the  energy  or  power  spectrum.  The  power  spectral  density  function 
specifies  the  contribution  of  the  harmonic  components  to  the  mean  square 
value  of  the  random  data  over  the  specified  frequency  range. 
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The  average  squared  value  of  the  random  data  will  approach  an 
exact  mean  square  value  as  the  observation  time,  T,  approaches  infinity 
for  a  frequency  range  between  f  and  A   f  (Ref.  2). 

In  equation  form: 

<%L£-ff+Af]=M^    ^//Z(ti-fAf)dt  _  (10) 

where  x(t,f,  A   f)  is  that  portion  of  x(t)  in  the  frequency  range 
bounded  by  f  and  f  +  A  f .   In  this  analysis,  of  course,  x(t)  could 
represent  either  O    or  w  .   For  small  values  of  Z\  f,  the  power 
spectral  density  may  be  computed  to  be: 

In  precise  form; 


(11) 


A-f-*0      *  *f->0  aT-o*MT/  >>  (12) 

o 

where  G  (f)  represents  the  power  spectral  density  function  of  either 
O      or  w  .   The  power  spectral  density  is  always  a  real  valued  non- 
negative  function. 

As  noted  earlier,  the  power  spectral  density  is  a  Fourier  trans- 
form of  the  autocorrelation  function  and  for  stationary  data  may  be 
expressed  as: 

Gx({)=zf/?x(y)ei2Tr  yJy=  4//?x(y)cos zttS  ydy 

-«/  o  (13) 

Thus,  from  Ref.  2,  the  mean  value  of  x(t)  is  given  by 
//  -//  <jy//C^/  »  an<*  tne  mean  8quare  value  is  given  by  0  =  /(^yCf)^' 
The  mean  square  value,  therefore,  is  equal  to  the  area  under  the  plot 
of  the  power  spectral  density  function  versus  frequency. 

From  Ref.  1,  typical  power  spectra  for  the  vertical  component 

of  turbulence  in  clear  air,  cumulus  clouds  and  thunderstorms  are 

2 
plotted  in  Fig.  11.  The  power  spectral  density  is  in  units  of  (ft/sec) 

per  radian/ft,  and  it  is  plotted  versus  the  reduced  frequency  in 

radians/ft,  with  an  additional  scale  of  wavelength  in  feet,  also  along 

the  abcissa.   Since  the  gust  disturbances  are  essentially  space 
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Figure  11.  Typical  power  spectra  of  vertical  component  of 

turbulence  measured  in  clear  air,  cumulus  clouds 
and  thunderstorm. 
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disturbances  rather  than  disturbances  in  time,  a  reduced  frequency, 
_/Z   ,  in  radians/ft,  is  normally  used  in  the  power  spectral  analysis 
of  turbulence.  The  relative  heights  of  the  curves  in  Fig.  11  indi- 
cate the  variation  of  the  intensity  of  the  turbulence  encountered 
under  the  different  weather  conditions.  The  similarity  of  the  slopes 
of  the  power  spectra  is  apparent,  and  their  magnitude  is  equal  to 

/-^   ,  where  C  is  a  constant.  The  spectra  show  a  characteristic 
rapid  decrease  in  power  with  Increasing  frequency. 

The  root  mean  square  value  of  the  gust  velocity  is  a  measure 
of  the  Intensity  of  the  turbulence  and  is  equal  to  the  area  under  the 
portion  of  the  spectrum  evaluated.   From  Fig.  11,  it  may  be  noted 
that  the  root  mean  square  value  of  the  turbulence  varies  in  magnitude 
from  4.48  ft/sec  in  clear  air  to  a  value  of  13.38  ft/sec  in  a  severe 
storm.  The  above  example  is  presented  to  describe  the  type  of  infor- 
mation that  Is  normally  available  from  the  power  spectral  density 
plot  of  a  given  parameter. 

The  power  spectral  density  of  O     and  w  are  plotted  on  pages 
80   and   85  of  Appendix  A.  As  expected,  the  power  spectral  density 
of  w  (series  2)  diminishes  with  increasing  frequency.  The  power 
spectral  density  plot  of  Q    (series  1)  also  decreases  with  increasing 
frequency;  however,  a  significant  mode  of  response  at  the  programmed 
short  period  frequency  of  .6  cycles  per  second  is  apparent. 

The  general  dependence  of  one  set  of  data  (such  as  &    )  upon 
another  (w  )  may  be  described  in  more  detail  using  the  cross  correla- 

O 

tion  function.  From  Ref.  2,  two  time  histories  of  data  are  illustrated 
in  Fig.  12. 


X(t)  A 


Figure  12.  -  Cross-correlation  measurement. 
(Reproduced  from  Ref.  2) 
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Like  the  autocorrelation  function,  the  average  product  of  x(t) 
at  time  t  and  y(t)  at  time  t  +  7  is  taken  over  the  observation  time 

a 

t.  For  this  analysis  x(t)  is  comparable  to  Q   ,  and  y(t)  is  comparable 
to  v  or  vice  versa.  As  T  becomes  infinite,  the  resulting  average 
product  approaches  an  exact  cross  correlation  function.  As  illus* ra- 
ted in  Ref.  2: 


/?xyCy)=  £ttt_^,/jCfd:yfttti/t 


r-*eo    '   /      u  (14) 

O 

R  (  "J   )  is  always  a  real  valued  function  which  may  be  either  posi- 
xy 

tive  or  negative;  however,  unlike  the  autocorrelation  for  random  data, 

/?.  C'y)        may  not  necessarily  have  a  maximum  at  7  *  0  or  be  an 

even  function  only.   R  (  y  )  displays  symmetry  about  the  ordinate 

xy 

when  x  and  y  are  interchanged,  so  that  R  (7  )"R  (-7)> 

xy       xy 

The  absolute  value  of  the  cross  correlation  function  is  bounded 
by: 

/R/y(y)/  <i/Rx(o)+Ry(o)/ 

where  R  is  comparable  to  the  autocorrelation  function  of  w  ,  R  is 

*  •        8  y 

comparable  to  the  autocorrelation  function  of  0    ,  and  R   is  the 

cross  correlation  function  based  on  Q     and  w  .   The  cross  correlation 

.8 
function  (labeled  as  cross  covariance)  for  &    and  w  is  tabulated  on 

g 
page  86   of  Appendix  A,  for  both  a  positive  1   and  a  negative  y   , 

and  plotted  versus  time  on  pages  87  and  88. 

The  Fourier  transform  of  the  cross  correlation  is  taken  to 
determine  the  cross  spectral  density  function  which  is  comparable  to 
the  power  spectral  density  function  for  a  single  variable.   Since  the 
cross  correlation  is  not  an  even  function,  the  cross  spectral  density 
function  is  a  complex  number  composed  of  the  co-spectral  density  func- 
tion (real  part)  and  the  quadrature  spectral  density  function  (imagin- 
ary part).  Therefore: 

Gxy  (f)  =  Cxy  (f)  -  /Qxy  (f)  <15> 
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where  C   Is  the  co-spectral  density  function  and  Q   is  the  quadra- 
ture spectral  density  function. 

The  co-spectral  density  function  is  the  average  product  of 
x(t)  and  y(t)  within  a  narrow  frequency  interval  between  f  and  f  + 
A.   f  divided  by  the  frequency  interval. 

The  quadrature  spectral  density  function  is  similar  to  the  co- 
spectral  density  function,  except  that  either  x(t)  or  y(t)  is  shifted 
in  time  to  produce  a  90  degree  phase  shift  at  the  frequency  f. 

The  co-spectral  density  function  is  expressed  as: 


o 
and  the  quadrature  spectral  density  function  is: 

7~ 

Q^^l^^^r^^A-ft/^^dt       '  (17) 

o 

where  x(t,f,^  f)  and  y(t,f,  A  f)  are  filtered  portions  of  x(t)  and 
y(t),  and  y*(t,f,  A  f)  signifies  a  90  degree  phase  shift  from  y(t,f,<Af) 
The  co- spectral  density  function  and  the  quadrature  spectral  density 
function  for  0    and  w  are  tabulated  on  page  89  of  Appendix  A. 

o 

The  cross  spectral  density  function  may  also  be  expressed  in 
complex  polar  form  as: 

GXy(f)=  /GXy(f)/<2  (18) 


The  magnitude/G  (f)/  of  the  cross  spectral  density  function  is 
equal  to: 

/Gx/  (f)l  =  VCxy  V)  *  Q\y  (f) '  <»> 

and  the  phase  of  the  cross  spectral  density  function  Q      (f),  is 

*y 

expressed  as:  #-  . 

6xy(f)=  fonJ_cxym\ 

(20) 
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The  magnitude  and  phase  of  the  cross  spectral  density  function 
are  tabulated  on  page  89  of  Appendix  A  and  plotted  on  pages  90  and 
91.   Again  a  significant  peak  in  the  amplitude  of  the  cross  spectral 
density  function  is  noted  in  the  region  of  the  programmed  short  period 
frequency  of  .6  cycles  per  second. 

The  amplitude  of  the  transfer  function,  between  the  input, 

x(t),  and  the  output,  y(t),  is  obtained  from  the  amplitude  of  the 

cross  spectral  density  function,  G   (f),  and  the  amplitude  of  the 

xy 

power  spectral  density  of  the  input,  G  (f).  The  amplitude  of  the 
transfer  function,  input  x(t)  to  output  y(t),  is  expressed  as  H  , 
and  is  equal  to: 


*        G. 


/-/v,   =    ^«yvL  (2i) 


In  this  analysis  x(t)  corresponds  to  w  and  y(t)  is  likewise  0      ', 
so  that: 

/Al=  ^lL  (22) 

The  phase  of  the  transfer  function  is  the  same  as  the  phase 

of  the  cross  spectral  density  function  as  outlined  above. 

& 
The  amplitude  and  phase  of  -r-  versus  frequency  are  tabulated 

on  page  92  of  Appendix  A,  where  the  phase  of  the  transfer  function  is 

expressed  as  portions  of  a  circle.   The  appropriate  translation  to 

angular  magnitude;  in  degrees,  is  noted  on  the  computer  program  beside 

the  printed  figures. 

Two  separate  curves  showing  the  phase  and  amplitude  of   — - 

a 
versus  frequency  were  plotted  from  the  results  of  the  power  spectral 

analysis  and  are  included  as  Figs.  13  and  14.   It  is  noted  that  the 

-TT-  .  flnd  also  a  phase  shift  of  -2_ 

>v„  *  ^ 

programmed  short  period  frequency  of  .6  cycles  per  second.   Although 
scatter  of  the  data  points  is  apparent  at  frequency  values  greater 
than  the  short  period  frequency,  the  amplitude  peak  of  the  transfer 
function  and  a  distinctive  phase  shift  may  be  readily  determined  from 
the  plotted  curves. 


peak  amplitude  of      ,  and  also  a  phase  shift  of  -2-  ,  occurs  at  the 
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frequency    /c/c/es   per  second) 
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THEORETICAL  SOLUTION  FOR  THE  PHASE  AND  AMPLITUDE  OF  THE  TRANSFER 

FUNCTION 


.  In  order  to  determine  the  accuracy  of  the  phase  and  amplitude 

of  tj-   ,  as  computed  by  the  power  spectrum  method,  the  appropriate 

longitudinal  equations  of  motion  for  the  F-100C  airplane  were  solved 

on  a  theoretical  basis.   Using  the  programmed  values  of  C    and 

mot 

Qn  ,  as  well  as  the  standard  stability  derivatives  for  the  F-100C; 
the  phase  and  amplitude  of  —  were  computed  over  the  same  frequency 
range  as  that  used  for  the  power  spectrum  analysis. 

The  following  basic  assumptions  were  made  in  determining  the 
equations  of  motion  and  solving  for  —  : 

1.  The  equations  of  motion  were  written  with  respect  to  the 
stability  axes. 

2.  The  aircraft's  elevator  was  considered  to  be  locked,  so 
that  there  were  only  three  degrees  of  freedom.   The  equa- 
tions were  written  therefore  for  linear  motion  along  the 
X  and  Z  axes  and  for  pitching  moments  about  the  Y  axis. 

3.  The  aircraft  was  considered  to  have  a  linear  response  to 
the  turbulence  encountered. 

k.     No  lateral/longitudinal  coupling  occurred. 

5.  The  airspeed  was  stabilized  at  the  beginning  of  the  run 
and  the  engine  thrust  was  constant  throughout  the  run. 

6.  Quasi-steady  flow  was  assumed. 

7.  Only  the  vertical  component  of  turbulence,  w  ,  contributed 

to  the  aircraft's  normal  acceleration. 

The  phase  and  amplitude  of  the  transfer  function  between  w 
•  g 

and  S     were  initially  obtained  using  a  simplified  solution  in  which 

the  phugoid  motion  was  considered  to  have  a  negligible  effect  upon  the 
short  period  oscillations  due  to  the  turbulence.   See  Ref.  9.   This 
simplified  method  is  designated  as  method  A  in  this  section.   To  check 
the  theoretical  results,  a  second  solution  involving  all  of  the  essen- 
tial terms  included  in  the  stability  equations  is  included.   This  exact 
method  is  designated  as  method  B.   Both  data  sheets  for  methods  A  and 
B  are  included  in  this  analysis  together  with  the  appropriate  values  of 
the  stability  derivatives  and  the  important  terms  used  in  each  solution. 
The  results  of  both  theoretical  solutions  were  in  close  agreement. 
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Q 

The  phase  and  amplitude  of  rj-   obtained  by  methods  A  and  B  are 
plotted  in  Figs.  15  and   16  .   It  can  be  noted  that  the  peak  amplitude 
and  the  phase  shift  of  -~      occurred  at  the  programmed  short  period  fre- 
quency of  .6  cycles  per  second. 

From  Section  13  of  Ref.  9,  the  longitudinal  equations  of  motion 
are  developed  with  respect  to  the  stability  axes: 


(23) 


(24) 


Xax/s:  u  *  ge  cos  ea  =  tu(ccs  f)u  +  T1m^6rps,cos$ 

2-  ax  '  s :  tv  -  C/0  v  r  g&  S/n  0o  -  -  Tu  (sin  f)cf-  T^^r^  S/n  $ 

+  ?uU +?%?  +  ?»»  +  ?*&  +*s£te  *-?srif  +*&& 
%+ch,r*.    Q^^2  Tuu+^m  -q^firn+MuU  +  Mtf 

These  equations  may  be  further  simplified  since  the  angle  0^     is 
assumed  to  be  small  so  that  sin  Qa   *  0  and  cos  G0    *  1.   In  addition, 
since  there  is  no  elevator  deflection,  flap  deflection,  or  speed 
brake  deflection  during  the  run,  S£  .  S r    and  Sg     all  equal  zero. 
Also  under  the  initial  conditions,  the  thrust  is  assumed  to  be  con- 
stant; therefore,  T  is  equal  to  zero,  and  the  variation  of  the  thrust 

with  RFM,  T  t„™,,  is  also  zero.   X  ,  the  variation  of  the  force  along 
3  RFM  q 

the  X  axis,  and  Z  ,  the  variation  of  the  force  along  the  Z  axis,  due  to 

q 

the  angular  pitching  velocity,  are  negligible.   In  addition  X  and  Z  , 

w      w 

the  variation  in  the  forces  along  the  X  and  Z  axes  due  to  the  time 
rate  of  change  of  the  vertical  velocity,  are  negligible.   For  this 
analysis,  the  vertical  velocity  is  considered  to  be  composed  of  two 
components,  w  and  w  ,  where  w  is  the  steady  state  component  of  the 

O 

vertical  velocity  and  w  is  the  fluctuation  of  the  vertical  velocity 
due  to  turbulence. 

Eliminating  the  terms  as  noted  above,  separating  the  vertical 
velocity  into  the  steady  state  and  turbulent  components,  and  intro- 
ducing La-Placiar  notation,  the  equations  of  motion  become: 
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su+ge    =  xuU  +  X^w  +  X^w9 


&  -u0se  =  2tJu  +  2u,ti  +  2„ Wj 


■u 


fv 


Re-arranging  the  equations  and  solving  for  the  w  terms: 

8 

(s-Xu)u        -xwu     +  ge     =  X^Wg 


(26) 
_(27) 
(28) 


(26) 


-(27) 


(28) 


for 


Using  Cramer's  rule,  the  above  determinant  is  modified  to  solve 
so  that  the  transfer  function  is  equal  to: 


(s-Xu)  -Xh  X* 

-2u        (S-  2  J)  &H 


(s-Xu)  -Xv*  f 


Expanding  the  determinant: 

_£     -       XHf?«(M&S  +  M„)  +  MH(s-g„$-2„£(s-Xu)(MlS+M») 
A/o  7 


(29a) 


-MuiJ]  +  M„f(s- Xu)(S-?w)  -  2u Xj] 


(29b) 
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where  A       denotes  the  original  determinant  in  the  denominator  of 
Equation  (29a).   Equation  (29b)  can  be  further  simplified  to: 

3  7-  ^  (M„  -XuM^)  -  Mw  (Xu  +  ?»)]s  <29c> 

A. 
Multiplying  the  numerator  of  the  ratio  of  the  determinants  by  8,  an 

expression  for  -S-  is  obtained,  so  that: 

/•o 

A    =    (Z^Mts  +  M^s3 '  +  fX*,2«.M*  ^XwMtu 

+  *v  (My  -XuM^)  -  Mw  (Xu  /  Z^IS2-  <30*> 

Method  A:   The  Simplified  Solution  of  the  Longitudinal  Equations 
of  Motion 

For  the  simplified  solution,  the  original  determinant  of  the 

denominator  is  considered  to  be  equivalent  to  the  product  of  the 

quadratic  expressions  for  the  short  period  and  phugoid  motions,  in 

Laplacian.  form,  so  that: 


£  =  feh,M&  +  M„)sJ+/'X»  ?u  M&+X„ Mu  +?„  R/  -Xu  Mw)-Mw  f/u  +Z„j]s 


"?  (s'-+Z?„ul„S+  ^Zf>«)(s^2J_,^S/>S+  lu*sp) 


(31) 

where  fPH    and  J^p     are  the  damping  ratios  of  the  phugoid  and  short 
period  motions  respectively.   And  6Jp#    and  CJsp      are  equal  to  2  7T 
times  the  phugoid  and  short  period  motion  frequencies.   The  contribu- 
tion of  (jdptf   is  considered  to  be  of  negligible  significance  in  the  air- 
craft's dynamic  response  to  turbulence  so  that  the  terms  containing 
UJp-H    can  be  eliminated.   Therefore,  the  denominator  can  be  approxi- 
mated by: 


(Sz)tiz  +  2  ?v£>OJSp  S  +  "J) 


7SP  ^sp  ^   f  w^ 

2 
Canceling  s   in  both  numerator  and  denominator  and  replacing  the 

operator  with  jw  ,  the  simplified  expression  for  —    becomes: 

i  -  (Zw  M^  +  Ml  v  )jiu  *■  Xw  ?u  M&  +  Xw  Mu.  +  Zw  (M„  -  Xu  M»)  -M  „  (Xu.  +  2 ») 
"3  (oS-sp—uj-*-)  +  -j2?spuJsf>u)  (32) 
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<v 


The  linear  velocity  along  the  X  axis  is  considered  to  be  constant  so 
that  the  perturbations  Z  ,  M  and  X   may  be  neglected.   Equation  (31) 
therefore  reduces  to: 


(33) 


The  amplitude  of  the  transfer  function  (method  A)  is: 


G 


(34) 


^(uJ3p-ou^+  4j/p  U>s%  uj2-' 
and   the  phase  of  the  transfer   function  is: 


*o_    =      90°-  tan"'    Z^^pUJspto 
**  fol-ou*)  (35) 

The  expressions   for  the  stability  derivatives  are  given  by 
Ref.    9  as: 


*,*    --^(^oc+C*) 


(36) 


"*  '  Iff  (Cn'*> 

Also   from  Ref.    10: 

CUsp    -    -LocMy-Moc  (39a) 

2  ?sp  CUsp  =l*~M%-Mk  ( AOa) 

where : 

^  * ,//  -  m  G«  e^  «*  ~*  ^  e(= v  <*-($ (42) 
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and   therefore: 

A^  =  flu.  €f,  ~/MiJ-  [2  fsp  CUspJ  (40b) 

AU=  -fi«e,jMf-£cos%]  <39b) 

The  computed  values  of  the  stability  derivatives  noted  above, 
as  well  as  the  quantities  used  to  compute  them,  are  listed  in  Table 
II.   Using  the  values  of  the  stability  derivatives  as  computed  in 

Table  II,  the  phase  and  amplitude  of  7—  are  computed  for  a  range  of 

h/g 

frequencies  from  .05  cycles  per  second  to  2.00  cycles  per  second.  The 
values  of  the  phase  and  amplitude  of  -^7-  ,  as  computed  by  the  simpli- 
fied solution  (method  A),  are  listed  in  Table  III.   The  phase  and 
amplitude  of  -—-    ,  as  computed  by  the  simplified  solution,  are  plotted 
versus  frequency  in  Figs.  15  and   16  . 

Method  B:   The  Exact  Solution  of  the  Longitudinal  Equations 
of  Motion. 

The  initial  computational  steps  of  the  simplified  solution, 

method  A,  as  well  as  the  same  given  stability  equations  and  initial 

assumptions  are  implied  for  the  exact  solution  as  well;  however,  for 

the  exact  solution  the  entire  determinant  for  w  will  be  retained  in 

f? 

the  denominator  of  the  ^-  term  throughout  the  analysis.   Starting 
with  Equation  (30a),—  is  equal  to: 

where  ^   is  equal  to  the  expansion  of  the  determinant  formed  by 
the  longitudinal  stability  equations,  solved  for  w  . 

A  =  (s-Xu)[s(S-Zh)(s-M<j)-SU»(m&S  +rtj)] 
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TABLE  II:   FLIGHT  DATA  AND  COMPUTED  STABILITY  DERIVATIVES 
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-sX„  Bu  (S-rt$+<}Mu  (s-2„)  -SUoX„  Mu  <43b> 

A  =  ^-(Xu+Z^  +  My+UoM^^+CMyfXu  +  Z^  +  Xui^-Xwitu. 

-  U«  (M„  -Xu  M^S2-*  [Uc  (XUM„  -X»  M„)  +  M%(Xu?u  -Xt,  ?w) 
+  q(M»+?«Mb}Js  +gUuM*-MuZJ)  (43c) 

For  simplification,  Che  real  and  imaginary  parts  of  the  numera- 
tor and  denominator  are  considered  separately  and  designated  as: 

R^     The  real  portion  of  the  denominator 

I      The  imaginary  portion  of  the  denominator 

R^     The  real  portion  of  the  numerator 

I      The  imaginary  portion  of  the  numerator 

N  o 

Q 
Therefore,  the  amplitude  of  the  transfer  function  - —  ,  is 

expressed  as: 


O  /    2.  2         ' 


ty      V  R^  +  i** 

and  the  phase  of  the  transfer  function  is: 


f0z  -f-ory-'Jji   -  fan''  *R.  <45> 
where: 

AW  =   '  z£n„(Xu+Z„)-Xw2uMi  -X^Mu  -2»(M»-XuMfrj]  (46) 

RN   •■    to2  A  (46a) 

7  N  -    UJ  7  ^  M  w  +  M„  7  (47  ) 

TN  U  *  £>  (47a) 

+  Mc(Xw2>>  -X:i'iv)  +  q(Mm-2uM{,)}>  (48) 

Jr        Ouit^l     *■  <--)  <*8a> 
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-  Uo  (M»  -  Xu  At  *)7  *   ?^  ^  ~  MtA  2w)  (49) 

/?^«    cu*  -  cozC    +- D  <49*> 

The  stability  derivatives  as  listed  for  method  A,  and  their 
values  as  tabulated  in  Table  II,  apply  to  method  B  as  well.   In  addi- 
tion, the  following  derivatives,  computed  in  accordance  with  Ref.  9, 
are  used  to  solve  for  the  coefficients  A  through  F,  of  Equations  (46a) 
through  (49a) : 

Xu   -  -f£*    CD  (SO),  *«   -  -g*   CL  (52) 

x»=  fgfo-cti  <*\     ^=-(-Brhc^  <53> 

The  computed  values  of  these  stability  derivatives,  as  well  as 
the  quantities  used  to  compute  them,  are  listed  in  Table  II.   The  values 
of  A,  B,  C,  D,  E,  and  F  are  listed  in  Table  IV.   The  phase  and  amplitude 

Q 

of  -—   are  computed  and  listed  in  Table  IV  for  a  frequency  range  from 
.05  cycles  per  second  to  2.00  cycles  per  second.  The  phase  and  ampli- 
tude  of  ~—   ,  as  computed  by  the  exact  solution,  are  also  plotted  versus 

fa/Q 

frequency,  together  with  the  values  obtained  from  method  A,  on  Figs.  15 

and  16.   The  two  different  theoretical  methods  produced  almost  identi- 

A 
cal  curves.   For  both  methods  the  peak  amplitude  of  ■—-  ,  as  well  as  the 

phase  shift,  occurred  at  the  programmed  short  period  frequency  of  0.6 

cycles  per  second. 
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DISCUSSION  OF  RESULTS 

Comparing  the  power  spectral  solution  and  the  theoretical  solu- 
tion (Figs.  13,  14,  15,  and  16)  by  matching  the  phase  and  amplitude 

A 
curves  for  -rr   ,  it  is  noted  that  in  both  cases  the  phase  shift  of  the 

transfer  function  is  in  close  agreement  and  the  greatest  amplitude  of 

Q 

— — •  ,  as  established  by  both  solutions,  occurs  at  the  same  frequency, 
0.6  cycles  per  second,  however  a  discrepancy  in  the  value  of  the  peak 
amplitude  of  -Z-  is  apparent.   The  amplitude  of  the  transfer  function, 

rvdl 

as  determined  by  the  power  spectral  method,  is  consistently  less,  by 
an  almost  constant  difference,  than  the  amplitude  as  established  by  the 
theoretical  solution  for  the  same  range  of  frequencies.  The  validity 
of  the  transfer  function  amplitude  as  established  theoretically  is  con- 
sidered to  be  accurate  inasmuch  as  two  different  solutions  produced 
identical  results  and  the  peak  amplitude  of  0.3  is  a  realistic  value 
for  an  aircraft  such  as  the  F-100.   Refer  to  Myslett  and  Smith  (10). 
The  power  spectral  solution,  as  derived,  would  imply  that  in  flight 
the  aircraft  was  not  as  disturbed  by  the  turbulence  as  the  theoreti- 
cal solution  predicted.  This  is  more  obvious  when  examining  the  auto- 
correlation  function  for  &    as  computed  and  graphed  by  the  digital 
program.   Although  the  distribution  of  the  fluctuations  of  the  auto- 
correlation  function  for  Q    was  in  agreement  with  theory;  the  magni- 
tude of  the  fluctuations  was  extremely  small  with  an  average  absolute 
value  of  .000001. 

The  computational  results  and  the  observation  mentioned  above 
led  to  a  close  scrutiny  of  the  instrument  sensitivity  on  the  test  air- 
craft used  to  gather  the  data.   In  particular,  the  sensitivities  of 
the  instruments  measuring  pitch  angle,  pitch  rate,  normal  acceleration 
at  the  C.G.,  and  angle  of  attack  were  checked  because  these  parameters 
provide  the  most  vital  information  regarding  the  actual  dynamic  response 
of  the  aircraft  to  turbulence.   The  entire  flight  record  for  Run  number 
8  was  re-examined  for  the  typical  inflight  measurements  of  the  para- 
meters mentioned  above.  The  results  are  summarized  in  Table  V. 
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Table  V:   Instrument  Sensitivity  and  In-Fllght  Measurements  for 
Fit.  26,  Run  8. 


Parameter  Measured 

Inst.  Sensitivity 
(In  terms  of  trace 
deflection) 

Average  Deflection 
Flight  Record 

Typical 

Resultant 

Value 

a 

12°/inch 

.22  inch 

tl.6° 

o 

e 

9.1° 
sec/ inch 

.07  inch 

+      1°/ 

.7  /sec 

Qt 

4g's/lnch 

.09  inch 

-  .36g 

<*•„ 

A. 5°/ inch 

.15  inch 
.05  inch 

.9*;<+^) 

.2°(-*w) 

The  tabulations  above  indicate  a  generally  small  response  of  the 
aircraft  to  the  clear  air  turbulence  encountered.  The  data  from 
Table  V,  the  diminished  magnitude  of  the  autocorrelation  function  for 

&      ,   and  the  consistent  trend  of  the  power  spectral  solution  for  the 

*  • 

amplitude  of  -~-  to  be  less  than  the  theoretical  amplitude  of  — - 
indicate  that  the  sensitivity  of  the  instrumentation  aboard  the  test 
aircraft  was  insufficient.  This  conclusion  was  further  verified  during 
the  investigation  of  the  instrument  sensitivity  at  NASA  when  discussion 
with  the  NASA  engineers  revealed  that  the  aircraft  used  in  collecting 
this  data  had  been  originally  instrumented  for  the  purpose  of  collecting 
maneuvering  type  data.  The  instrument  sensitivity  requirements  for  this 
type  of  maneuvering  flight  are  not  as  stringent  as  those  necessary  for 
the  power  spectral  type  analysis  of  turbulence.  The  accuracy  of  the 
response  is  even  more  critical  in  a  power  spectral  analysis  since  the 
associated  digital  computation  is  of  such  length  that  any  initial  dis- 
crepancies will  be  magnified  by  "rounding  off"  errors  which  are  carried 
throughout  the  computer  program.  The  amplitude  discrepancy  of  the 
transfer  function  noted  above,  however,  was  of  indirect  value  in  an 
academic  sense  in  that  it  generated  a  thorough  secondary  analysis  of  the 
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flight  test  data  and  a  better  appreciation  of  the  entire  flight  test 
program  from  considerations  of  the  instrument  location  and  sensitivity, 
to  the  digital  program  itself. 
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SUMMARY  AND  CONCLUSIONS 

Despite  the  single  discrepancy  in  the  solution  for  the  transfer 
function,  — —  ,  it  is  felt  that  this  analysis  did  represent  a  signi- 
flcant  and  unique  application  of  power  spectral  methods  because  the 
transfer  function  between  a  random  input  and  the  resultant  aircraft 
dynamic  response  was  established  on  a  continuous  basis  over  a  finite 
frequency  range.  The  amplitude  discrepancy  between  the  theoretical 
and  power  spectral  solutions  for  the  transfer  function  is  considered 
to  be  entirely  due  to  improper  instrument  sensitivity  and  not  to  the 
mathematical  computations  or  the  power  spectral  analysis  itself.  A 
closer  correlation  between  the  amplitudes  of  the  transfer  functions 
would  be  realized  using  the  same  analytical  procedures  if  the  in- 
flight response  measurements  were  more  accurate. 

Future  effort  in  turbulence  response  analysis,  particularly  with 
power  spectral  applications,  will  require  that  more  sensitive  instru- 
mentation be  incorporated  on  the  test  aircraft.   In  addition,  the  angle 
of  attack  vane  should  be  re-located  from  its  present  position  on  the 
vertical  stabilizer  to  a  position  at  least  two  and  a  half  fuselage 
diameters  ahead  of  the  aircraft  as  proposed  in  Refs.  1  and  11. 

Although  this  analysis  considered  a  singular  application  of 
the  power  spectral  method,  additional  areas  of  investigation  related  to 
gust  loading  would  be  worthwhile  for  future  research.  A  preliminary 
study,  by  Press  and  Houbolt  (12),  has  been  conducted  to  determine  the 
variation  in  normal  acceleration  with  center  of  gravity  position  for  a 
fighter  type  aircraft.  Through  power  spectral  analysis,  the  ratio  of 
the  root  mean  square  values  of  the  normal  acceleration,  determined  by 
the  area  under  the  gust  spectrum,  was  readily  established  between  center 
of  gravity  positions.   It  was  determined  that  as  the  center  of  gravity 
position  is  shifted  aft,  the  peak  frequency  response  moves  to  a  lower 
frequency  and  diminishes  slightly.   Other  significant  areas  of  future 
study  and  investigation  include  the  effect  of  low  short  period  damping 
on  the  gust  loads  of  aircraft  in  continuous  turbulence  and  the  magni- 
fication in  gust  load  stresses  that  results  from  wing  flexibility. 
Refer  to  Press  and  Tukey  (13). 
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In  general,  the  power  spectral  methods  of  analysis  may  be  re- 
garded as  an  effective  and  accurate  tool  for  analysis  of  a  random  type 
disturbance.   In  particular,  it  is  a  powerful  method  for  further  ex- 
ploration of  the  science  of  gust  loading. 
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Appendix  A 


The  Sequence  of  Digital  Computer  Steps  in  the  Power  Spectrum  Analysis 

The  mathematical  steps  outlined  in  Table  VI  explain  the  compu- 
tational procedure  behind  each  of  the  digital  computer  steps  in  the 
power  spectrum  analysis  (Ref.  14).   All  information  in  Table  VI  refers 
to  Fig.  17,  and  each  step  in  the  table  is  related  to  the  appropriate 
computer  operation  by  the  computer  sequence  number. 

The  symbols  used  in  Table  VI  are  defined  as  follows: 

X  i   value  of  discrete,  equi-spaced  time  series  X(t) 

X  i   value  of  time  series  X(t)  after  subtracting  mean 

Y  i   value  of  time  series  Y(t)  after  subtracting  mean 

P  Lag 

n  Number  of  discrete  data  points 

m  Maximum  lag 

t  Constant  time  interval 

R^P'   Autocovariance  of  series  X  at  lag  p,  R^P'  -  R  *P-At' 

X  °  r*        X        X 

P      Power  spectral  estimate  of  series  X  at  frequency 

x  m  a  t 

„(h)    _  h  it 

P    -  P  vrzre 
x       x 


SP     Smoothed  power  spectral  estimate  of  series  X  at  frequency 

h  V 


xy 


*y 


m  a  t 


(p) 
R  r    Cross-covariance  between  series  X  and  Y  at  lag  p 
xy 

P      Cross-power  spectral  estimate  at  frequency  h  77 


m  A  t 


C      Cospectrum,  the  real  part  of  P^  '  at  frequency  h  7r 
Xy  xy  m  A   t 

Cy     Quadrature  spectrum,  the  imaginary  part  of  P^  '  at  frequency 
xy  xy 

h  7T 


m  A  t 


SC  Smoothed   cospectrum  at   frequency  h   n 


m  A   t 
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SO/  '  Smoothed  quadrature  spectrum  at  frequency  h  ^ 

xy  mAt 

AM  Amplitude  of  cross  spectrum  at  frequency  h  7f 

xy  m  *  t 

PHAS  Phase  of  cross -spectrum  at  frequency  h  7f 

xy  m  ^  t 

T  Transfer  function  at  frequency  h  rr 

xy  m  ^  t 

TAM  Amplitude  of  transfer  function  at  frequency  h  7f 

xy  m  A  t 
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TABLE  VI:   Sequence  of  Digital  Computer  Steps  In  The  Power  Spectrum 

Analysis 


COMPUTER 
SEQUENCE 
NUMBER 


MATHEMATICAL  OPERATION 


REFER  TO  PAGE: 


© 


The  mean  of  the  data  is  first 
calculated  and  subtracted  from  each 
value  of  X' . 


»>*=££.  *i 


'x      n 


-c-/ 


(Original  Data) 


e      76 


X*.  -  x*. "  ^x 


i=  /tzr..%n 


ho  81 


© 


The  autocovariance  is  then  computed. 

(p)        f        1-& 


e 

A/c 


7 


77,78 
82,83 


® 


The  raw  estimate  of  the  power  spectrum 
is  obtained  by 


777 


/      77       A-o 


O  < p<  m 


e 


A/c 


NOTE:^^^)^^) 


79 


84 


© 


The  raw  estimates  of  the  power  spectrum 
are  then  smoothed  by  "hamming" 


SPX«*  =  .54Px<o)  +  .46Px 


(h) 


(h-i) 


(o) 
(f>). 


0) 


e 


(h+l) 


5PX  =.Z3PX     +.54%   +.23% ;    o<h< 


tn 


Wn 


86^ 


85 
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TABLE  VI    (Continued) 


COMPUTER 
SEQUENCE 
NUMBER 

MATHEMATICAL  OPERATION 

REFER  TO  PAGE: 

© 

A  check  sum  is  computed  to  check  the 
computations  of  the  estimates. 

This  should  equal  R   ,  the  autocovari- 
ance  at  zero  lag. 

b      79 

Wj           84 

© 

The  cross -covariance  is  computed  by  the 
formula 

n-p 

s?(-p)-  '  Tv     y    D 

^«y   "^V?'  P=o,i,z...,m 

86,87, 
88 

■ 

© 

©   ! 
© 

The  cross-spectrum, /P   ,  is  given  by 

/^y     =  Cxy      +  x.  Oxy 

The  cospectrum  is  obtained  from 

The  quadrature  spectrum  is  obtained  from 

QxV   =  n  Z,  €p(Ryy  +R/y  J5/n-£r 

[J          O  <  p<m 

Both  the  cospectrum  and  quadrature 
spectrum  are  smoothed  by  "hamming"  as 
shown  in  step  4. 

89 
89 

73 


TABLE  VI  (Continued) 


COMPUTER 
SEQUENCE 
NUMBER 


MATHEMATICAL  OPERATION 


REFER  TO  PAGE: 


The  amplitude  and  phase  of  the  cross- 
spectrum  are  computed. 


a^ix^p  [5o(:r 


89,90 


© 


(h) 


(h) 

OQxy 


RHA5v    --  ARCTANi 


The  phase  angle  computed  is  the  phase 
shift  of  series  Y  with  respect  to 
series  X. 


89,91 


© 


The  transfer  function  from  series  X 
to  Y  is  given  by 

(h)  ,a   .      (h) 


Txy     ~  px(h) 


TAM/y  C 


92 


© 


The  amplitude  is  given  by 


(h) 

rAMxy  =  AtAxv 


The  phase  shift  from  X  to  Y  is  the  same 
as  that  of  the  cross-spectrum.  Both 
the  transfer  function  amplitudes  from 
X  to  Y  and  Y  to  X  are  printed. 


92 
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TMC02T-AUT0C0VAR UNCETHJB   POWER   5PECTR7KANALY5I5  -   WBIW  OF   MR.    ».    1W 
HEALTH    SCI6NCES    COMPUTING    FACILITY,    UCLA 


P  H~0   6   L    E   M      nUHIER 


RUN     I 


INPUT  DATA  TO  BE  PKINTEO  OUT YES 

/rj.,ei*Jr      ic,     •Wc-wd 

INPUT  SERIES  TO  BE  PL0TT60  OUT ' 

DE TRENDING 

PREWHITENING 

VALUE  ST  CONSTANT  CT^Ed  TN  THE  FKEBIIITEIINB  TRANSFORMATION  Z(T)«X<T*I>-Cm)   -—     -0. 

NUMBBR  OF  SERIES  -   2 

NUMBER  OF  DATA  POINTS  •   920 


NUMBER  OF  LAGS  CHOSEN  -  100 
NUMBER  OF  SELECTION  CARDS  •   1 


USE  PR6VI0US  DATA 


CONSTANT  TIME  INTERVAL 0.10000  SECOND 

NUMBER  OF  VARIABLE  FORMAT  CARDS  •    1 


-& 


ORIGINAL    DAT*   OF    SERIES 


0.O0O27 

0.00052 

-0.00023 


0.00027 

0.00102 

-0.00023 


0.00140 
-0.00023 
-0.00023 

0.00127 
-0.00023 
-0.00073 

0.00421 

0.00152 

-o.oooii 

0.00516 
-O.0OS*9 

-0.0003* 
0.00729 
0.00090 
0.0006* 


0.00077 
0. 00115 

-0.00023 


0.00027 
0.00002 
0.00202 

-0.00023 
0.0001* 
0. 00115 

-C.00C23 


0.00002 

0.00290 

-0.00023 


0.00027 
0.00002 

0.00315 
0.000** 


0.00002 
0.003*0 
0.00190 


0.00002 
0.00152 
0.00313 


0.00002 
0.0006* 

0.00*03 


0.0002 J 


O.OOOOt 
O.OOOOt 
0.00371 


0.00027 

0.00002 

-0.00023 

0.00278 


-0.00023 

0.00115 

-0.00099 


-0.00111 

0.00190 

-0.00136 


-0.00211 

0.00177 

-0.002*9 


-0.00211 
0.00090 
-0.00g*9 


0.00265 
-0.00023 
-0.00066 


0. 

,00278 

0 

.00152 

-0, 

.00161 

0 

.00667 

-0, 

.0012* 

-0, 

,00036 

0.008*2 
0.00090 

o.aons 


0.00027 

0.00027 

-0.00023 

0.0C01* 

-0.00237 

-0.001U 

-0.00086 

-0.00086 

-0.00099 

-0.00099 

-As-OOJl* 

-0.00*87 

0.00202 

-O.OOC23 

-0.00161 

0.00278 

0.00152 

-<U003*9_ 

0.00667 

C.0001* 

-0.0003* 

0.00830 

O.C0090 

0.00115 

0.00027 

0.0006* 

-0.00061 

-0.00086 

-0.00099 

-0. 00312. 

-0.00061 

-0.00036 

-C.O.0il*_ 

-0.00211 

-0.00099 

0.00039 

-0.002*9 

-0.00199 

-0.00036 


0.00152 

-0. 00023 
-0.002J1 


0.00077 
-0.08023 
-0.00199 


0.00215 

0.00340 

-0.00387 


0.001*0 
0.00*16 
-0.00317 


0.00039 
0.00052 

-0. 00199 
0.00039 
0.00414 

-0.00142 


-0.00023 
0.00002 

-0.00161 
0.00039 
0.0065* 

-0.0003* 


-0.00211 

0.00002 

-0.00249 


-0.00149 

-0.00014 
-0.00040 


-0.00048 
-0.00011 

0.0004* 


-0.00023 
-0.00099 
-0.00136 


-0.00024 
-O.OOOIA 
-O.OOOU 


-0.00023 

-0.00023 

0.00240 


0.00039 
0.00634 
0.00039 


0.00064 
0.00420 
0.00177 


0.00115 
0.00165 
0.00403 


0. 00604 
0.00O1* 
-0.00036 


0.00233 
0.00014 
-0.O0034 


-0.000*1 

0.00014 

-0.00036 


-0.00211 
0.00014 
-0.00034 


-0.00437 
-0.00061 

-0.00036 


-0.00*33 

-o.ooo6a 

C. 00077 


-0.00*37 

-0.00036 

0.00*16 


0.00840 
0.00090 

0.00102 


0.00027 
0.00039 
-JUOOOil 
-0.00086 
-0.00099 
-0. 003.12. 
-0.0012* 

-o.coon 
_=o_,OQa?j! 
-0.0012* 
-0.00099 
0.00039 
-0.00106 
-0.00186 
-0.001*9 


0.00729 
0.0006* 

0.00052 


0.00*03 
0.0001* 
0.00027 


0.00152 
0.0001* 
0.00027 


0.0001* 
0.0001* 
0.00027 


-0.00016 
0.0001* 
0.00027 


0.00002 
0.0001* 
0.00027 


-0.00073 
-0.00036 
-0.00011 


-0.00086 
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FREQUENCY                 (J)  POWER  SPECTRAL  ESTIMATES 
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FREQUENCY 

CO-SPECTRUM 

QUA0RA1URE  SPECTRUM 

AMPLITUDE  OF  CROSS-SPECTRUM 

PHASE  OF  croSs-SpCCTruM 

ICYCLES/SECONOI 

OF 

OF 

OF 

OF 

SERIES   1  AND   2 

SERIES   1  AND  2 

SERIES   1  AND   2 

SERIES   1  ANO   2  • 

-0. 

-0.0021232 

-0.0032647 

0.0039111 

0.6566696 

0.050 

-0.0014642 
-0.0007291 

-0.0036339 
-0.0018874 

0.0039178 

0.6690397 

0.100 

0.0020233 

0.6913315 

0.150 

-0.0002058 

-0.0016613 

0.0016740 

0.7303843 

0.200 

0.0000841 
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-0.0015001 

0.0015024 

0.7589148 

0.250 

-0.0013435 

0.0013496 
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0.300 
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-0.0009608 

0.0009824 
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0.500 
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0.0002839 
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0.0008557 
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0.1393655 
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0.1451075 

1.050 
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0.0001651 
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0.0000451 
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"     Tr.mamsz 
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Q^\      GRAPH  OF  AMPLITUOE  OF  CROSS-SPECTRUM  OF  SERIES   1  AND   2  (IN  LOG  SCALE)  AGAINST  FREQUENCY  ICYCLES/SECONO I 
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(J]}   PHASE  OF  CROSS-SPECTRUH  Of  SERIES   1  AND   2—  PLOTTED  IN  FRACTIONS  OF  *  CIRCLE  AGAINST  FAEOUENCY  (CYCLES/SECONO) 
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(S) 

FREQUENCY 

COHERENCE  SQUARE 

AMPLITUDE  OF 

AMPLITUDE  OF 

PHASE  OF 

ftlASe    or 

(CYCLES/SECOND) 

OF 

TRANSFER   FUNCTION 

TRANSFER   FUNCTION 

TRANSFER   FUNCTION 

T-BAMs/'**  rv*<:r'o*J 

SERIES   1  ANO   2 

JRJJH   I  TQ_  1 

FROM   2  TO   1 

FRON   1  TO   2  • 

rAo"*    f   To  Z    (aboi*e*) 

-0. 

0.0619165 

6074.7177734 

0.0000102 

0.6586698 

IZ3  O 

0.0500 

0.1614117 
0.5300088 

7526.7287598 

0.0000214 

0.6890397 

//  Z.I 

0.1000 

4609.3800659 

0.0001  ISO 

0.6913315 

/ll  .1 

0.1500 

0.6889446 

3483.0234375 

0.0001978 
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